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The DFT/B3LYP/6-311G(d,p) study of water and methanol hydrogen 
bonded clusters is presented. The effect of hydrogen bond formation on 
the geometry of clusters is considered. With increasing cluster size 
shortening of the hydrogen bond as large as = 0.2 À is observed. The 
accurate binding energies, corrected for the basis set superposition error, 
are reported. They reflect strong cooperative effects accompanying the 
hydrogen bond formation. The maximal value of binding energy per 
hydrogen bond is observed for cyclic hexamers. Comparison with earlier 
studies on water clusters is carried out. The obtained results are discussed 
in terms of the cluster populations in the bulk phase. 


1. INTRODUCTION 


The hydrogen bonded liquids, such as water, alcohols, etc., reveal specific 
spectroscopic and thermodynamical properties. They are commonly described 
in terms of the thermal equilibrium of hydrogen bonded molecular clusters of 
different sizes. For example, elongation of the O-H bond, accompanying the 
hydrogen bond formation, causes the red shift of the O-H bond stretch 
frequency, usually by more than 100 cm”. The upfield migration of the 
hydroxyl proton NMR signal in alcohols with temperature is attributed to 
breaking of the hydrogen bonds — small clusters with strongly shielded 
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hydroxyl protons are favoured at high temperature. Also the enhanced 
association of molecules (monomers) at lower temperatures, causes the 
structure of a liquid to be more compact — the density of a liquid increases. 

Proper description of physical features of hydrogen bonded liquids requires 
detailed knowledge of the properties of individual clusters. They can be 
computed with high accuracy with modern ab initio quantum mechanical 
programs. Water clusters were already subjected to theoretical investigations; 
a number of them were published over the last five years alone. For example, 
the structure and energetics of (H;O),, n=12,16,20 clusters was examined by 
molecular dynamics simulations [1]. Structure, energetics and IR spectra of 
medium-sized clusters, (H2O),, mn=7-10, were also computed [2-4]; 
the calculation included ab initio MP2 method in the determination of 
molecular geometry. The MP2 structures, binding energies and IR spectra of 
cage (H2O),, n=11,12 clusters were published recently [5]. The calculations 
were carried out with double zeta plus polarization basis set and included 
corrections for the basis set superposition error. Weinhold [6] considered 
aseries of water clusters and reported their properties: electronic energies 
corrected for the basis set superposition error (BSSE), rotational and vibrational 
temperatures and zero point energy (ZPE) corrections. These data were used in 
a framework of the Quantum Cluster Equilibrium (QCE) theory [7] to compute 
the populations of clusters in liquid water as a function of temperature. Similar 
calculations were carried out for other hydrogen bonded systems, including 
N-methylacetamide [8], ethanol [9,10], etc. These authors were able to 
reproduce qualitatively a number of thermodynamical and spectroscopic 
features (including variation of the density, heat capacity, chemical shifts, 
quadrupole coupling constants, etc. with temperature). However, the low level 
of theory used to compute the clusters’ properties (i.e. Self Consistent Field 
(SCF) method and small basis set) is one of the factors responsible for only 
qualitative agreement between experimental and theoretical data. 

In this paper we present the detailed Density Functional Theory (DFT) [11] 
study of molecular clusters of water and methanol. Clusters (H20)p, 
n=1,2,3,4,5,6,7,8,12 and (CH3OH),, n=1,2,3,4,5,6,7 are considered. The limited 
set of clusters was chosen on account of their expected dominance in the liquid 
phase over the range of temperatures [6,9]. The calculations include: molecular 
geometry optimization, zero point energy corrections and basis set 
superposition error corrections. The average energy per hydrogen bond for all 
clusters are also reported and discussed in terms of the cooperativity of 
hydrogen bonding interactions. The presented results can be used in order to 
predict the composition of liquid water and methanol within the QCE model. 
We believe, that the determination of cluster populations from more accurate 
data will provide a new insight in the structure of the hydrogen bonded liquids. 
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It will be also helpful in the rigorous judgement on the predictive capabilities of 
the QCE model. We will address these problems in the forthcoming papers. In 
Section 2 the details of DFT calculations are reported. In Sections 3 and 4 the 
results for isolated monomers and clusters are described and discussed. In 
Section 5 we state our conclusions. 


2. COMPUTATIONAL PROCEDURE 


It is well known that the simplest ab initio method — the SCF method — often 
fails to predict accurately the equilibrium geometry of molecules in 
near-degenerate ground states. On the other hand reliable geometries for typical 
closed-shell systems can be obtained already at the SCF level of theory, 
provided the relevant basis set is chosen. In our previous paper [12] we 
investigated ethanol molecule — the well-behaved, closed-shell system. 
However, the geometry of the C-O-H fragment turned out to be somewhat 
inaccurate at the SCF level. Correlated ab initio methods (e.g. MCSCF, MP2, 
etc.) are usually superior to SCF calculations. Accurate molecular 
conformations are also predicted by the Density Functional Theory [13], in 
particular when the B3LYP exchange-correlation potential [11,14,15] is used. 
For the ethanol molecule [12] — the system similar to those considered 
here — the method reproduces well the experimental (microwave) data. In this 
paper molecular geometries were optimized at the DFT/B3LYP level using the 
6-311G(d,p) basis set [16] for all water and methanol clusters. The zero point 
energy corrections were computed within the harmonic approximation by the 
standard FG—matrix method [17] for all clusters at the same level of theory. 

The binding energies of the M, clusters ( M =H,0, CH;OH), 


AE, — nE non - E, (1) 


were corrected for the well known basis set superposition error (overestimating 
binding energies in small basis sets) using counterpoise (CP) correction 
method [18]: the energy of each monomer u in the cluster M, is computed 


(i) using the basis set of the whole cluster — Ej", and (ii) using the basis set of 


the monomer alone — Ef”. The CP correction, defined as 


6° = YD ES”) (2) 
is positive and includes the effect of a small distortion of a monomer in the 
cluster. The counterpoise corrected binding energy is thus given by the 
following formula: 
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AES = AE, -— 6° (3) 
The percentage cooperative enhancement (relative to dimer) was computed as 
ae a 
(% coop), = AES (4) 


where m, denotes the number of the hydrogen bonds in the cluster. 


All calculations reported in this paper were performed with the parallel 
version of the ab initio quantum chemistry PQS package [19,20] implemented 
for a cluster of Pentium II/400 MHz PC machines. The maximum speedup 
obtained with 12 CPUs was about 11 as compared to the single processor run. 


3. THE WATER AND METHANOL MONOMERS 


In Tab. 1 theoretical and experimental equilibrium geometries of water and 
methanol molecules are compared. Bond lengths and valence angles predicted 
by the DFT/B3LYP/6-311G(d,p) method are very close to the microwave 
spectroscopy data. The relative error usually amounts to less than 1%, 
indicating high accuracy of the reported values. It should be emphasized, that 
the molecular geometry (as well as other properties: rotational and vibrational 
temperatures, etc.) used in the QCE modeling of liquid water structure by 
Weinhold [6], were obtained at much lower level of theory (i.e. SCF/3-21G). 
Consequently, the bond lengths were too short (rotational temperatures 
inaccurate) and vibrational frequencies (and temperatures) too high. We may 
also expect substantial errors in Weinhold's binding energies of hydrogen 
bonded clusters. As we shall see they are overstabilized at this level by as much 
as 20-30%. It is clear, that the determination of reliable cluster populations 
should be based on much more accurate than SCF/3-21G data. 


Tab. 1. Theoretical and Item geometries of water and methanol molecules 


Bond m] (in À) Valence [eee (in 9) 


Water Theory Expt.^ 0. — EE] — __ 
Bran RSEN 9575 104.51 


[| CO ^ CH’ OH | C-O-H H-C-H 
mid Expt.“ | 1 Tm 1.0906 0.961 107.84 108.27 
1.425 1.094 0.945 108.53 108.63 


“See eg. [24] 
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4. THE WATER AND METHANOL CLUSTERS 


The ab initio based QCE calculations of populations in liquid water [7] 
included 9 water clusters, the largest of which was octamer. This model was 
extended: 18-cluster QCE treatment [21] comprised ring, tetrahedral and 
fullerene-like clusters up to 26-mers. Nevertheless all the calculations of 
clusters' properties were carried out at the simple SCF/3-21G level of theory. 
For this reason only qualitative significance should be given to the predicted 
populations of water clusters in the bulk phase. 

In this paper we report results of a series of DFT/B3LYP calculations on 
water clusters (up to 12—mers) and methanol clusters (up to 7—mers). Geometry 
optimization and harmonic vibrational frequencies calculations (at the 
optimized geometries) were carried out at the DFT/B3LYP/6-31 1G(d,p) level. 
All frequencies turned out to be real indicating that considered structures are 
stable (they correspond to local minima at the potential energy surface). 
Energetics, including ZPE and BSSE corrections, and the most important 
geometrical parameters are presented in Tables 2 and 3. The structures are 
shown in Figures 1 and 2. 

Geometry optimization usually led to the cyclic clusters. Only in the case of 
water and methanol dimers, linear structure is more stable than the 
corresponding cyclic structure; the latter is in fact transition state due to high 
strain in the four-membered ring. Only one linear cluster — the methanol 
pentamer — was found to be stable. However, its total binding energy, Eq. (1), is 
by about 9 kcal/mol higher than for the cyclic counterpart. It is hard to judge at 
this point whether the global minima on the potential energy surfaces were 
found. It should be emphasized that the geometry optimization of hydrogen 
bonded species is a difficult task due to the presence of numerous local minima 
and transition states. Our experience with similar calculations on ethanol 
clusters [12] leads to the conclusion that local minima of cyclic structures for 
a given number of monomers in the cluster are very close in energy. In fact they 
have the same topology: the oxygen-hydrogen ring is preserved and the only 
geometrical differences are in the relative orientation of the hydrofobic C-C 
chains. As discussed elsewhere [7] it is sufficient to include only one 
representative $n$—mer in the QCE calculations. This is why we did not carry 
out a search for global minima of water and methanol clusters. 

In all water and methanol clusters the lenghts of the free OH bonds (rg); 


refer to Figures 1b' and 2b' for details), i.e. bonds that do not participate in 
hydrogen bonding, remain essentially the same as they were in corresponding 
monomers (= 0.96 À). 


J. Jaroniec, K. Ziólkowska, P. Borowski, T. Janowski and K. Woliñski 


190 


'(.t [ 910311) si9gurexəu/siourenə) jo seue[d ət UseMI0q SZU puoq 93e1oAV , 

'(.1 [ AmI) sig9gurexəu/srəurenə1 jo soue[d oq ur su18uə[ puoq IFLIAY acy 

p) ‘bap 

`191sn|o eq ur spuoq uoSoJpÁq Jo 1oquinN. 
"pouoda: are v rid puoq o8e1oA? ‘(q | Ə1m3r,T) sjuourge1 

O-H: : -O 3} ur puoq HO pue puoq uasospAy jo 1red e jo sy3uo| — 804 pue [Hu ‘puog HO 221J ? Jo qiaue[ — Rog 
`uonoərtoəo Ág1ouo 

urod 013Z — 9dZ '459 — “TV = 454 V (C) ‘by "10119 uontsodrədns 19$ siseq — ,59 (T) ‘by “(8rəuə Surpurq —" V, 


90€ v6T '9| S6ç'SIT| 8T€'SG | 8276 | %¿I9'0¿I Ieureoopoqq 


Lev 621 š '9| T0£'c4 | £98'¿€ | 6£T'6 | F99'601 
ISO ITI j "4| 96V'ES | S¿€'¿¿ | ESS TT | T2808 
608°€6 j L9V' LE | IG87T | 988'0I | 8I€'S9 
S¿T 8 ; LI@LE | ¿GT9T | 699'01 | HHE'ES 
89€'c9 i "4| GOL8%@ | SLO €I | PEC OL | SLT Tv 
9cv'9y G| O@6'ST | I89'8 | 0028 | 009'F2 
580 97 i v69TvI | c8€8 | 6972 |926'€6 [SUL L 
62683 j 99/' I48'6 | LE9L |2£97 ruq 
08€'€I Jourouo]q 


(y ut sun8uəl[ puoq *,[our/meoy ut sors1ouo) (d'P)O I L€-9/d X TE 8/.LG :s1o1snpo 103^ AL `Z ARL 


I 


NUVO 
1eure)doH 
IourexoH 
Jourgjuoq 
Ioureno], 
c-Tou, 


6 
¿ 
9 
S 
y 
€ 
€ 


191 


Structure and energetics of the hydrogen-bonded water and methanol clusters 


`(p) bap 

“1ə1sn[o IP ut spuoq uoSoJpÁuy Jo 19quinN.; 
“poyodal are sulšuə[ puoq o8esaae ‘spuog QD pue HD » sqduo[ — 094 pue H2, *(.q Z Ə1mn8r1) 

suauge O—H: - -O IP ur puoq HO pue puoq uo3o1pÁu JO 1red e jo sui3uo[ — E04 pue H ( 94 *puoq HO 221 ? Jo ui3uo[ — i 
"UOrjoalioo 

K81ouo utod 0132 - JdZ ‘459 — “AV = aS V '(C) “bg 10119 uontsodsadns yas siseq — 439 XD "bg *K81ouo Surpurq -" y, 


q 


689'v£c 9IZ'IS | 9€E°0% | £6c'0T| cS02Z | L 1oure1dəH 
98€ 102 COV Ph | E88 LI | I6£'OT| OPE | 9 C-1ourexoH 
876 TOZ £9c'£v |9€€'8I | £0£'OT | 618 19 | 9 [T-39urexəH 
979291 LG9°S€ | 6L9 FI | L29001) SECOS | G C-19urjuoq 
160791 65926 | IPL EL | 8P£'Ot| O6E IP | v [-1ourejuod 
9I8'€£I TEG LZ | GIT IL | G8€6 | ChESe | F Joule nel 
0ZT00T crcl |0£€ | 1697 | PS6%Z] € JOWUL 
196°0 | 6%2°¢9 PLSP (SPL | 9TLL | 914 I nuq 
196°0 | S60°ZE IIWOUON 


pia u 


(y ut sujSuo[ puoq *,[ourjeox ur səršrəuə) (d'p)ODI [£-9/4X'1C8/1L4d :siə)snio jouemjow au], `ç AQL 


192  J.Jaroniec, K. Ziólkowska, P. Borowski, T. Janowski and K. Woliñski 


Fig. 1. Structures of water clusters: (a) monomer, (b) and (b') dimer, (c) trimer-1, (d) 
trimer-2, (e) tetramer, (f) pentamer, (g) hexamer, (heptamer, (i) and (i') octamer, (j) 
dodecamer 


Also the CO and CH bond lengths (<ç and ro) in methanol do not vary with 


cluster size: they are equal to = 1.42 À and = 1.095 À respectively. 
The association of monomers to form larger clusters results in profound 
changes in the geometry of O"H-O fragments. With increasing cluster size the 
length of a hydrogen bond ( r2) ) decreases — from = 1.92 A (dimer) to = 1.70 A 
(higher cyclic clusters) for water and from = 1.88 À (dimer) to = 1.69 À 
(heptamer) for methanol. This shortening of a hydrogen bond by about = 0.2 À 
is accompanied by elongation (by about = 0.02 À) of the OH bond in the 
O^H-O fragment. It should be noted, that incorporation of a large number of 
water molecules in the cluster results in spatial rather than planar structures 
(see eg. Figures li and 1j), containing hydrogen bonds of different lengths. 
Indeed, the iterative geometry optimization of water octamer, with initial 
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geometry of cyclic type, led to the cubic structure shown in Figure 1 (i). Note in 
passing, that cyclic octamer was incorrectly predicted to be stable at the 
SCF/3-21G level of theory [6]. The described structural trends are 
characteristic to all hydrogen bondend systems [6,9,12]. 


(a) 


Fig. 2. Strucutres of methanol clusters: (a) monomer, (b) and (b') dimer, (c) trimer, 
(d) tetramer, (e) pentamer-1, (f) pentamer-2, (g) hexamer-1, (h) hexamer-2, (i) heptamer 
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The content of a particular cluster in the bulk phase is closely related to its 
binding energy (AE) per hydrogen bond, i.e. AE/m,. The so-called 
percentage cooperative enhancement relative to dimer (% coop, Eq. 4), is 
a convenient measure of binding energy. Note, that even with basis set as large 
as 6-311G(d,p), the accurate determination of AE requires estimation of basis 
set superposition error. The BSSE correction for the systems considered in this 
paper constitutes 30-40% of AE. For this reason the counterpoise corrected 
[18] DFT/B3LYP/6-311G(d,p) binding energies (AE , Eq. 3) were used in the 
calculations of % coop. It is evident from Tables 2 and 3 that the binding 
energy per hydrogen bond (as well as % coop) exhibits the maximum value for 
medium-sized water and methanol clusters, i.e. cyclic hexamers in both cases. 
Such a behaviour is due to the nonpairwise—additive (cooperative) character of 
hydrogen bonding [22]. As described elsewhere (see eg. [23]) the hydrogen 
bond formation in dimer increases both the polarity of OH bond and basicity of 
the proton acceptor (the second water or methanol molecule) favouring 
formation of another hydrogen bond — the trimer structure is obtained. In the 
case of water and methanol clusters the cooperative effect is large — the binding 
energy per hydrogen bond of cyclic hexamers is by more than 60% higher than 
in dimers. 

The cooperative effect in water calculated with DFT/B3LYP/6-311G(d,p) 
method is somewhat lower than that computed earlier [6] at the SCF/3-21G 
level of theory. Also DFT binding energies are by 1—3 kcal/mol lower than 
those predicted by SCF. Such discrepancies may influence significantly the 
populations of clusters in the bulk phase. It can be easily understood if we 
recall that they are governed by the canonical partition functions of clusters [7] 


which are proportional to 
CP 
q” ss exp AE, 
ji kT 


It is obvious, that even small changes in the binding energy result in large 
changes in q“; this may significantly alter the cluster populations in the bulk 


phase. At the first sight pentamers, hexamers and heptamers should dominate at 
room temperature in liquid water and methanol (cf. Tables 2 and 3). It should 
be realized, however, that the described enthalpic effect is not the only effect 
responsible for the composition of hydrogen bonded liquids [6,7]. Numerous 
accessible rotational states and low energy vibrational states enlarge the 
entropy of a system. This may efficiently exclude some of the clusters in favour 
of others. The rigorous treatment of both enthalpic and entropic effects is 
provided by the QCE theory [7]. Implementation of QCE code is one of our 
main goals in the nearest future. 
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5. CONCLUSIONS 


We calculated a number of properties of the water and methanol clusters 
using the DFT/B3LYP/6-311G(d,p) approach. The calculations included: 
molecular geometry optimization, binding energies corrected for the basis set 
superposition error and zero point energy corrections. 

The bond lengths of water and methanol monomers are very close to the 
experimental (microwave) values. With increasing cluster size we observe 
substantial shortening (by about 0.2 A) of the hydrogen bond. This is 
accompanied by elongation of the OH bond in the O^H-O fragment by about 
0.02 À. The remaining geometrical parameters of clusters: the lengths of a free 
OH bond (in water and methanol), CO and CH bonds (in methanol), as well as 
some of the valence angles are essentially the same as in the corresponding 
monomers. 

The calculated DFT/B3LYP/6-311G(d,p) binding energies per hydrogen 
bond, ranging from 4.77 kcal/mol (water) and $4.57 kcal/mol (methanol) for 
dimer to 7.91 kcal/mol (water) and 7.41 kcal/mol (methanol) for cyclic 
hexamer, indicate strong cooperative effects in both liquids. They are lower 
than these obtained at the SCF/3-21G level of theory by about 1-3 kcal/mol. We 
expect, that this new energetics is more accurate than reported in the QCE 
treatment of liquid water and will provide more accurate cluster populations in 
the bulk phase. The problem of evaluation of the composition of hydrogen 
bonded liquids will be undertaken in our future papers. 
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